import tables
import pysam
import sys
import getopt
import numpy as np
import os
import dspomic.ptpileup

def errlog(*msgl):
  msg = ' '.join([str(x) for x in msgl])
  sys.stderr.write(msg + os.linesep)

def main():
  # some variables
  samFiles = [None]
  h5fname = "pileup.h5"
  downsample = 1.0
  mapQuality = 0
  removeDuplicates = False

  # get command line options
  try:
    opts, args = getopt.getopt(sys.argv[1:], "o:d:rq:")
  except getopt.GetoptError, err:
    errlog(err)
    sys.exit(2)
  
  for o, a in opts:
    if o == "-o":
      h5fname = a
    elif o == "-d":
      try:
        downsample = float(a)
      except ValueError, err:
        errlog(err)
        sys.exit(2)
      if downsample > 1 or downsample < 0:
        downsample = 1.0
    elif o == "-r":
      removeDuplicates = True
    elif o == "-q":
      try:
        mapQuality = int(a)
      except ValueError, err:
        errlog(err)
        sys.exit(2)
      if mapQuality < 0:
        mapQuality = 0
  samFiles = args
  h5ob = dspomic.ptpileup.PTPileup(h5fname)
  for fname in samFiles:
    h5ob.importIntervalFile(fname, ds = downsample, mq = mapQuality, rd = removeDuplicates, fmt = "bam")
  h5ob._close()
  
if __name__ == "__main__":
  main()